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Abstract 

It is known that the dependence of the emission-Une luminosity of a typical cloud 
in the active galactic nuclei (AGN) broad-line regions (BLRs)upon the incident flux of 
ionizing continuum can be nonlinear. We study how this nonlinearity can be taken into 
account in estimating the size of the BLR by means of the "reverberation" methods. 
We show that the BLR size estimates obtained by cross-correlation of emission-line 
and continuum light curves can be much (up to an order of magnitude) less than 
the values obtained by reverberation modelling. This is demonstrated by means of 
numerical cross-correlation and reverberation experiments with model continuum flares 
and emission-line transfer functions and by means of practical reverberation modelling 
of the observed optical spectral variability of NGC 4151. The time behaviour of NGC 
4151 in the Ha and H/3 lines is modelled on the basis of the observational data by 
Kaspi et al (1996) and the theoretical BLR model by Shevchenko (1984, 1985). The 
values of the BLR parameters are estimated that allow to judge on the size and physical 
characteristics of the BLR. The small size of the BLR, as determined by the cross- 
correlation method from the data of Kaspi et al. (1996), is shown to be an artifact of 
this method. So, the hypothesis that the BLR size varies in time is not necessitated 
by the observational data. 

Key words: galaxies: active — galaxies: nuclei — galaxies: Seyfert - galaxies: individual: 
NGC 4151 



1 Introduction 



In the early 1970's, in the course of observations of rapid variability of the optical spec- 
trum of the Seyfert galaxy NGC 4151, the time lag of va riations in the Ha line w i th re- 



spect to variations in the optica l continuum was discovered (iLvutvi fc Cherepashchukl. Il971 



Cherepashchuk fc Lyutyil . Il973l ) . The time lag was interpreted by iLyutyil (119771 . Il982l ) as a 
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consequence of the fact tha t the emission-line clouds a re at some distance from the ionizing 
radiation source. Later on, lAntonucci fc CohenI (119831 ) observed much smaller time lags in 
variations of the H/? and H7 lines. Such a difference in the time lag, in Ha greater than in 
is observed i n other active galac t ic nucl ei (AGN) as well (see Table 6 in lPeterson et al.ll2004l ). 
According to lShevchenkd (jl984l . ll985al ). this difference in time lag is due to an essential non- 
linearity in the dependence of the Ha luminosity of an individual cloud upon the ionizing 
continuum flux incident on the cloud, the dependence in the higher order Balmer lines being 
close to linear. This explanation was made with the assumption that the duration of the 
emission-line flare is much greater than the duration of the flare in the ionizing continuum, 
and the duration of the latter one allows its description by (5-function. However, continuum 
variations are not so rapid usually; therefore, in order to extract physical information from 
the observed emission-line variations, it is necessary, in addition to theoretical estimates, to 
use numerical modelling taking into account the timescale of continuum variations. 

In the present paper, we study theoretically how the nonlinearity in the emission-line 
luminosity, Li, of the broad-line region (BLR) cloud, in its dependence on the ionizing con- 
tinuum flux, Fi, incident on the cloud, can be taken into account in estimating the BLR sizes 
by means of the "reverberation" methods. We show that the BLR size estimates obtained 
by straightforward cross-correlation of emission-line and continuum light curves can be much 
(up to an order of magnitude) less than those obtained by reverberation modelling. First of 
all, we demonstrate this by means of abstract representative numerical cross-correlation and 
reverberation experiments with model continuum flares and emission-line transfer functions. 
Then we accomplish practical numerical modelling of the light curves of NGC 4151 in the Ha 
and H/5 lines on the basis the observational da ta by iKaspi et al.l (119961 ) and the theoretical 
model of the BLR by IShevchenkd (11984 . Il985al ). This model is characterized by allowing for 
thick geometries of the BLR, taking into account the anisotropy of line emission of individual 
clouds and, most important, taking into account the nonlinearity of the "Li-Fi" relation. 
This nonlinearity allows one to explain the differences in the time lags for different lines. 
Cross-correlation estimates of the BLR size are also made. They turn out to be small in 
comparison to the estimates obtained by the direct reverberation modelling. 

The values of the parameters of the BLR model are derived directly from the rever- 
beration modelling, and that is why we do not use any specific numeric results of modern 
photoionization models of the emission-line spectra of AGN. We use only the fact that, ac- 
cording to these models, the emission-line response of an individual cloud (particularly, in 
Ha and B.(3) can be nonlinear. We also allow for a constant emission-line component. 

Though the presented theoretical inferences on the time lags can be of general interest, 
our primary goal is to apply them to explaining the time behaviour of the Ha and H/5 
lines in the emission-line spectrum of NGC 4151, in the framework of a simple uniform 
one-component model. Of course, explaining the time behaviour of the whole emission-line 
spectrum would require much more complicated multi-component models. 



2 Effective Stratification of a Homogeneous BLR 



The broad-line r egion (BLR) o f an active galactic nucleus, according to the "standard 
model", see, e.g. |Petersonlll988l . represents an aggregate of line-emitting clouds under the 
effect of ionizing radiation of the central source. 

The dependence of the emission-line luminosity L\ of an individual cloud upon the value 
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of the incident ionizing flux Fi, in accordance with tfie pliotoionization m odels of spectra of 



Kwanlll984l : iMushotzky fc Ferlandlll984l ). is described by a 



active galactic nuclei (see, e.g. 
power law: L\ oc , where s > 0. 

The rate of heat input in a gas cloud optically thick in the ionizi ng con t inuum is directly 
proportional to the value of the ionizing flux incident on this cloud. iKwan (Il984l ) noted that 
therefore the cloud's emission-line luminosity should be, in a first approximation, directly 
proportional to the ionizing flux; different lines, however, behave differently. For example, in 
the case of La the dependence is somewhat weaker than linear. The La quanta leaving the 
cloud are produced in the traditional H ll zone. At the high ionization parameters typical of 
AGNs the collisional ionization from the excited levels of hydrogen (in particular, from the 
second level) are effective even in the Hll zone. With increasing ionization parameter their 
efficiency grows, a nd this leads to weakening of the specified dependence ( IKwan fc Krolikl . 



19811 : iKwanl . 119841 ). 



So, the emission-line response of an individual cloud can be nonlinear. This fact was 
recognized already in the first successful photoionization models. According to them, quanta 
in many lines are produced mainly not in the traditional H ll zone, but deeper, in the so-called 
"deep partly ionized zone" . Successful modelling of stationary optical emis sion-line spectra 
of AG N requires the following two circumstances to be taken into account (IKwan fc Krolikl . 



19791 ): the power-law shape of the spectrum of ionizing continuum (i.e., the fact that the 
major fraction of ionizing quanta is in the X-ray part of spectrum) and the big column 
densities of the clouds emitting in lines. If the X-ray luminosity of the ionizing source is 
great enough in comparison with the UV one, the "deep partly ionized zone" is formed in 
the cloud. Taking into account the contribution of this zone increases the luminosity of the 
cloud in Balmer lines, whereas the luminosity in La is stabilized at the level of the luminosity 
of the H II zone. So, the collisional amplification of Balmer and Paschen lines takes place in 
the "deep partly ionized zone". Inside this zone the excitation temperatures of these lines 
increase with optical depth, but ultimately attain some limiting values. The limiting values 
are insensitive to variation of the ionizing flux , because the Balmer and Paschen co ntinua 
dominate in cooling at such depths (IKwanl . Il984j ) . According to iKwan fc Krolikl (Il98ll ) , when 
collisional ionization becomes the main source of ionizations and cooling, the rate of cooling 
increases with increasi ng electron temperatur e approximately as exp(— 32 x 10'*K/Te). In 
the standard model by iKwan &: Krolikl (1l98ll ). Tg ~ 8000 K in "the deep zone"; the steep 
dependence of the rate of cooling on temperature, as Kwan and Krolik noted, provides only 
weak variation of Tg with depth and insensitivity of Tg to variation of the model parameters, 
in particular, the ionization parameter. Increasing the ionizing flux makes higher levels of 
hydrogen attain the limiting excitation temperatures; the luminosity of the "deep zone" in 
the relevant lines then ceases to react to changes of the ionizing flux, i.e., in this limit they 
are constant. In the Balmer series, the approach to the limiting temperatures affects first 
of all the Ha line, th e n H/3, and so on. Thus, according to the photoionization models 
(IKwan fc Krolikl . Il98ll : iKwanl . Il984l ). the dependence of the cloud's emission-line luminosity 
on the incident ionizing flux for the Ha line is weaker than for H/5, for H/5 is weaker than 
for H7, and so on. 

Due to the difference between Ha and H/? in the value of the s parameter, the Balmer 
decrement increases with increasing distance of clouds away from the central source, therefore 
a photographic BLR image (if such an image could be obtained) would be larger in Ha 
than in H/3. A formula for an effective BLR radius in a line wi th an arbitrary s v alue in 
the homogeneous model of the cloud aggregate was deduced in (jShevchenkol . Il985al ). This 
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effective stratification is explained by difference s betw een emission lines in the degree of 
nonlinearity of the Li{Fi) function. IShevchenkd (119881 ) showed that within the framework 
of the homogeneous model of the cloud aggregate, if one takes into account the re sults of 
the photoionization calculati o ns of the emission-line spectra of active galactic nuclei (IKwanl . 



iMushotzkv fc Ferlandl . Essl), it is possible to explain the observed time lags and 



amplitudes of variations in major optical and ultraviolet (UV) emission lines in the spectrum 
of NGC 4151. 

After the first successes of the photoionization computations of the AGN emission -line 
spectra, significant progre ss was made in this field; see, e.g., reviews by iFerlandl (120031) and 



Leighlv fc Casebeerl (|2007r). Multi-component models were proposed and studied (jCoUin-Souffrin et al. 



19881 : CoUin-Souffrin &: Lasota . 1988 : Korista et al. . 1997 ). which allowed to reproduce the 
relative fluxes in high-ionization and low-ionization li nes simultaneously. Evidence was 
found for the presence of opt i cally thin line-emitting gas ( Ferland. Korista &: Peterson . 1990l : 



Shields. Ferland &: Petersonl . Il995l ). This progress promoted much deeper understanding of 
the AGN emission-line spectra — it turned out that the uniform models are too simple to 
reproduce the whole spectra. However, in what follows, our study concerns only Balmer 
lines. We aim to explain the time behaviour of the Ha and H/3 lines in the emission-line 
spectrum of NGC 4151, in the framework of a simple uniform one-component model. Of 
course, explaining the behaviour of the emission-line spectrum in total may require much 
more complicated mult i- component models. 

Effective, not physical, stratification is present in our one-component model, due to the 
nonlinearity in each cloud's line emission. The alternative to a homogeneous BLR with 
effective stratificatio n is a physi c ally s tratified BLR. Investigating variability of the UV 
lines of NGC 4151, luirich et al.l (Il984[ ) offered the BLR model consisting of three zones 
with different physical characteristics (see Table 2 in their article). iGaskell fc Sparkd (119861 ) 
proposed a model consisting of two zones (see Table 1 in their article). These models are 
not considered henceforth; we adopt the effective stratification picture as imphed by the 
nonlinearity in cloud's line emission. 



3 The Reverberation Model 



Blandford fc McKed (119821 ) offered a procedure to recover the BLR structure by analysis of 
line and continuum lightcurves. This is the so-called method of "reverberation mapping" . Its 
essence consists in the following: the observed light curve in a hne is supposed to represent 
a convolution of two curves: the transfer function describing physical characteristics and the 
geometry of the BLR and the light curve in ionizing continuum. The emission-line luminosity 
of an individua l cloud was su pposed to depend linearly on the incident ionizing flux. 



Shevchenkd (119841 . Il985al ) found necessary and sufficient conditions for the existence of a 



time lag of a maximum of an emission-line flare in relation to a (short duration) continuum 
flare when the BLR structure is isotropic with respect to the central source; these conditions 
are: the typical cloud should emit in the line mainly from the side facing the central source, 
and, either a central cavity should be effectively present in BLR, or the s parameter in the 
formula L\ oc Ff should be less than one. These conditions set useful reference points for 
our modelling. In the 1990 's the effect of nonlinear response as well as of anisotropy of the 
individual cloud emissio n in specific BL R m odels, as applied to the cross-cor r elation anal ysis, 
was studied in detail by ISparkd (1l993l ) and lO'Brien. Goad fc Gondhalekarl (1l994j . Il995l ). 
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We adopt the homogeneous model of the cloud aggregate ( IShevchenkd . Il984l . Il985al ). 
The effective BLR radius R is defined by screening of the peripheral part of the aggregate by 
the clouds situated closer to its centre: R = {crn)~^, where a [cm^] is the mean geometrical 
cloud section orthogonal to direction to the central source, n [cm~^] is the cloud concentration 
(number of clouds in a unit volume). Generally, the BLR can contain a central cloud- free 
cavity of radius Rq. Let us remark that the Balmer quanta, unlike the ionizing quanta, can 
leave the BLR freely even at large cloud-covering factors of the "sky" of the central source, 
because the dispersion of the cloud velocities in the BLR is assumed to be great; the latter 
fact is testified by the large width of the observed emission lines. 

The model of a homogeneous (outside the central cavity) distribution of the clouds is 
equivalent, in what concerns the transfer function form, to a model with zero covering factor 
but with n oc e~^^^, an exponential decrease in the cloud concentration with increasing 
distance from the centre. In both interpretations, R characterizes the BLR radius for all 
lines. In the first case, it is the radius of the "lit" zone in the homogeneous aggregate, and 
in the second case it is the e-folding scale of the cloud concentration. 

We assume that a typical BLR cloud represents a flat "pancake" emitting lines solely from 
the side facin g the ionizing radiat ion source, and, what is more, emitting orthotropically. See 
discussion in (IShevchenkd . Il985bl ) on the physical basis for this assumption. The planes of the 
clouds are either orthogonal to the direction to the central source, or are oriented randomly. 
The phase function, describing the phase angle dependence of the cloud's line emission, is 
different in these two cases. For random orientation, the effective phase function (the phase 
function of a volume unit containing naany cl ouds) coincides with the phase function of a 
spherical cloud, provided (Shevchenko, 1985bh : the cloud is completely opaque in the line, 
the line quanta are produced at small optical depths, and the cloud surface emits in the line 
orthotropically. 

One should make a reservation that the pancake shaped cloud, as well as a uniform cloud 
aggregate model itself, is a physical idealization that can be used only as an approximation for 
the real arrangement of line-emitting material in t he BLR. The real structure might be closer 
to a combinat i on of a disk and an outflowing wind (lEmmering. Blandford fc Shlosmaru . Il992 ; 



Murrav et all . Il995l : IChiang fc Murravl . Il996l : iBottorff et all . Il997l : lElvid . [20001). We adopt 



the pancake shape fo r the BLR cloud exc lusively for convenience in mathematical modelling: 
indeed, according to IShevchenkd (1l985bl ). the phase function of the pancake cloud with the 
plane orthogonal to the ionizing source direction provides the maximum anisotropy of line 
emission, i.e., this is a physical limit worth theoretical examination, while the phase function 
of a spherical cloud (or, equivalently, randomly oriented pancakes) gives an approximation 
for the phase function of randomly oriented optically thick line-emitting material, i.e., it 
describes a situation that is expected to be closer to reahty. 

If the "pancakes" are orthogonal to the central source direction, the transfer function 
representing the dependence of the observed integrated emission-line flux f{t) on time t 
counted from the moment of the (5(t)-fiare of the central source in continuum, is as follows 
fIShevchenkd . Il984l . 11985^ ): 
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J g{r, t) dr, Ro<t< 2i?o, 
t 

Jg{r,t)dr, t > 2Ro, 



t/2 



where 



- 1 



^l-2Sg-r/R 



and r, R, Rq are measured in the hght-travel time units. 

In the case when the planes of clouds are oriented randomly, their mean phas e function co- 
incides with the phase function of a spherical cloud. This function is as follows (IShevchenkd . 
1985bh : 



j{9) oc (1 + cos( 



1 + - cos 6* 



(2) 



where 6 is the "ion izing source - clou d - observer" angle, 0<^^<7r,0<s<2. Gen- 
eral formula (3) in (jShevchenkd . Il984l ) for the transfer function, after substitution of phase 
function (12]), becomes: 



/Woe 



t 

R 



Jg{r,t)dr, 0<t<2Ro 



Ro 



(3) 



^Jg{r,t)dr, t>2Ro 



where 



t/2 



1 + 



t 



Transfer functions and can be expressed through incomplete 7 functions. The be- 
haviour of the transfer functions with different values of the s parameter (while Rq = 0, 
R = 15 It-days) is demonstrated in Fig. [H The qualitative difference in the behaviour of 
the functions with s less and greater than unity is clearly seen. In particular, f{t) peaks at 
t = for s > 1 and at t > for s < 1. 

The model emission-line light curve is determined by the convolution formula: 



roo 

F,{t) = a f{T)F:{t-r)dr 
Jo 



(4) 



where Fi is the integrated flux in the line, Fc is the observed flux in continuum, a is the 
normalizing dimensional factor. Since it is the optical continuum, not the ionizing one, that 
is observed, we use an assumption, formulated below, on a relation between the continua. A 
difference of expr ession (HI) from those usua lly u sed (valid in the c ase of the linear "Li-Fi" 
relation; see, e.g., Blandford fc McKee 1982 and Horne et al.l 2004 ) consists in raising of Fc 
to the power s. Let us remark that, according to (jll), on taking /(r) in the form of 5-function, 
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Figure 1: Transfer functions ([T]) and ([3]) (figures (a) and (b), respectively) for three values 
of the s parameter, while Rq = 0, R = 15 It-days. 
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one gets F](t) oc F ^(t), i.e., the dependence for the case of quasi-stationary spectrum; see 
flShevchenkol . ll988h . 



4 The time lag and the cross-correlation method 



Techniques for cross-correlation analysis of AGN emission-line variability have demonstrated 
remarkable progress during the last decade. The methods of calculation of the basic proper- 
ties of the cross-correlation function (CCF), namely, the lags of the CCF peak and CCF cen- 



1995, 


1998; 


Welsh 


1999) 



the continuum var iability time scale is important for correct estimation of the BLR size (e.g.. 



Edelson fc Krohklll988i l 



Let us consider the time lag as determined by means of cross-correlation analysis in the 
case of nonlinear emission-line response of an individual cloud. In this Section, we measure 
the time lag in model numerical experiments and study the dependence of the time lag on the 
parameters of a model transfer function and duration of the continuum flares. We consider 
the case of a single flare of various durations. As the model transfer function we take Eq. ([T]) 
corresponding to the case of the "pancake" clouds orthogonal to the central source direction. 
The central cavity in the BLR is set to be absent: Rq = 0- 

The model light curve in the continuum is assumed to have the form of the bell-like 
function Fc{t) = sech((t — to)/T), where to = 50 d and T is effective duration of the flare, t 
is time in days. The model emission-line light curves are computed on the time interval of 
500 d with the step of 0.05 d. 

In the upper part of Fig. El the model light curve in the continuum (T = 1 d) and the 
computed emission-line light curve are presented. The latter curve has been obtained by 
means of convolution of the light curve in the continuum and the transfer function ([1]) with 
R = 15 It-days, Rq = 0, s = 1. In the lower part of Fig. [21 the normalized cross-correlation 
function of these curves is plotted. The shift of the peak of the cross-correlation function 
is clearly visible; as determined numerically, Atpeak — 3.8 d. Note that the exponential-like 
decay of the resulting curves in both plots reflects the radial structure of the line-emitting 
region, and the rise in the CCF reflects the shape of the continuum flare. 

One may argue that the linear response model is just a linearization of the nonlinear 
response model, and that any BLR radius estimate made in the linear response model is 
therefore an approximation that might be not far from reality. However, one should take 
into account, firstly, that the "equal time-travel" paraboloidal surface inside the BLR covers 
a whole range of distances from the ionizing source just after the ionizing flare, secondly, that 
with increasing time after the flare this surface retreats from the source to larger distances. 
The slopes of the linearized dependences for individual clouds on the surface vary significantly 
in both space and time. The response slopes might be averaged on the surface, but then 
the change of the averaged slope with time should be taken into account; the latter is never 
done in practice. So, it is not surprising that the nonlinear response model, as compared 
to the linear one, can give very different quantitative results on the BLR radius. This 
directly follows from the qualitative differences in the response function for different values 
of the s parameter, as seen in Fig. [H A vivid manifestation of the insufficiency of the 
linearized response model is that increasing the time lag value in the nonlinear response 
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Figure 2: The upper plot: the model light curve in the continuum (the continuous line; 
R = 15 It-days, i?o = 0, T = 1 d, s = 1) and the computed emission-line light curve (the 
dashed line). The lower plot: the normalized cross-correlation function. 
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Figure 3: The dependence of the time lag on the parameter s; Rq = 0, T = 1 d. The straight 
hne segments represent theoretical relations ([5]). 



model can be achieved either by increasing the BLR radius or by specific increasing the 
response nonlinearity, namely, by decreasing the s parameter value in relation to unity (see 
rel ation (|^ ) , while in the lin ear response model only the BLR radius can be varied. 



Shevchenkd (jl985al . Il994l ) obtained an approximate theoretical relation of the time lag of 
the maximum of the emission-line light curve to the s parameter in the homogeneous model 
of the cloud aggregate with or without a central cavity (i?o > 0). The continuum fiare was 
described by a 5 function. This relation is as follows: 

_( Wil- s)R, 0<s<l- 2R,/{WR), 

\ 2Ro, s > l-2Ro/{WR), ^ ^ 

where the constant W depends on the choice of phase function; W = 3.19 in the considered 
case of clouds with regular orientation. In the case of phase function ([2]) one has W = 2. 

By examining the shifts of peaks of cross-correlation functions at various values of param- 
eters we can assess how well relation works when the ionizing fiare has a finite duration. 
Consider first how the time lag varies with s at fixed R (Fig. [3]). The time lag, as defined 
here, is the value of the distance along the time axis from t = up to the first maximum of 
the cross-correlation function; i.e., it is Atpeak- We do not examine the shift of the centroid 
of the cross-correlation function here. The parameter s is varied from 0.1 to 2.0 with the 
step of 0.01. The ionizing fiare duration is fixed at T = 1 d. The curves for i? = 5, 15 
and 25 It-days are plotted. Theoretical dependences for R = 5, 15 and 25 It-days and 
Rq = are plotted as straight line segments. It is clear that the theoretical and numerical 
estimates of the time lag at such a relatively small duration of the continuum fiare are in a 
good agreement. 
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Figure 4: The dependence of the time lag on the BLR radius R; Rq = 0, T = 1 d. 



Now we consider the dependence of the time lag on the BLR radius R (Fig. H]). The 
ionizing flare duration is the same, T = 1 d. The dependences for s = 0.5 and s = 1 are 
plotted. We see that the linear character of theoretical relation ([5]), valid for the S{t) ionizing 
flare, is preserved in the case of s = 0.5. In what concerns the s = 1 case, it is completely 
different. From the viewpoint of relation this case is degenerate, and the predicted 
value of Atpeak is constant (zero). In reality we observe the "Atpeak^-R" relation similar to a 
logarithmic one. 

So, as follows from Fig. [3] and Fig. HI for s > 1 the CCF Atpeak value depends only weakly 
on R. The cross-correlation peak time lag for s > 1 is small in comparison with R expressed 
in the light-travel time units. The dif ference can reach a n order of magnitude. For s = 1, a 
similar phenomenon was observed by IPerez et al.l (119921 ) for the CCF centroid estimates of 
the BLR size in an isotropic (with respect to the ionizing source) model of the line-emitting 
cloud distribution. They found that such estimates can be less than the real size of the BLR, 
the difference reaching two times. This is in accord with our findings. 

What is the role of the ionizing flare duration in the degenerate case s = 1? The 
dependences of the time lag on R in this case for various fixed values of T are plotted in 
Fig. [5l We see that this role is far greater than that of the BLR radius. The observed 
"Atpeak^-R" dependences seem to be logarithmic indeed; one can verify that, unlike the 
rational and power-law functions, the functions of the form a + b\n{R + c), where a, b, c are 
fitting parameters, provide ideal visual description of the observed curves. Note that further 
work is required to explore how the logarithmic dependence found here is vulnerable to the 
choice of the ionizing flare shape. 

In summary, according to our numerical-experimental findings, the BLR radius is only 
weakly related to the measured Atpeak value in the mathematically degenerate but observa- 
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Figure 5: The dependence of the time lag on the BLR radius R at various values of duration 
T of the ionizing flare; s = 1, Rq = 0. 



tionally most common case of s = 1. The role of the timescale of variability is far greater. 
Therefore, the lines with s ~ 1 are of little help in determining the BLR size by means of 
cross-correlation techniques; instead, the lines with s essentially less than 1, such as Ha, 
should be used for this purpose. This conclusion has been obtained in a model framework 
and thus may be model-dependent in some way. However, it makes clear that there are no 
general theoretical grounds to believe that Atpeak is mostly determined by the BLR size. 
Note that solely the CCF peak offset has been examined. The CCF centroid offset should 
be examined as well in a future study to check whether it exhibits the same behaviour. 

5 Reverberation modelling of the emission-line light 
curves of NGC 4151 

In this Section we examine the effect of taking into account the nonlinearity in a cloud's 
line emission in practical modelling of emission-line variability of an AGN. We model the 



emissi on- line variability of the Seyfert galaxy NGC 4151. The observational data of lKaspi et al 



(Il996l ) on variabilit y of the nucleus of t his galaxy in Ha, H/? and optical continuum is used. 



The observations of lKaspi et al.l (119961 ). performed in the framework the AGN Watch Pro- 
gramme, cover the time interval of approximately three months, from 1993 November until 
1994 February. While we use the light curve continuum data for the whole time span of the 
observations, the data on the fluxes in the lines during the first five days and during the last 
five days of the observatio nal time span are e xcluded, following the usual practice of eliminat- 



ing the border effects (see Maoz et al. ( 199ll )). The error bars of the individual observations. 
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defining the weights of the observations, are taken into account in the modelling. 

We optimize t f ie mo del parameters by using a non-linear least-squares method (iLevenbergl . 
19441 : iMarquardtl . llQGSh to minimize y ^, thereby finding best-fit parameter values and their 



standard errors (see Press et al. f 19971 )). We find that the iterations converge to the same 



solution for random starting values in ranges specified below, with the deepest minimum of 
in all cases, though several other local minima exist, as demonstrated below. We expand 
the standard errors by the square root of the reduced x^, because the best-fit reduced is 
greater than unity in all best-fit models found. 

As the transfer functions, we use Eqs. ([T]) and (jHl)- The best approximation of the 
computed model light curve to the observed one is found by means of the modelling. We 
vary six parameters: the BLR radius the nonlinearity parameters Sopt for Ha and H/3, 
the normalizing coefficients a for Ha and H/5, the flux F^n in the narrow component of Ha 
(on the F^n value for H/3 see below). 

The radius Rq of the central cavity (a zone free from line-emitting clouds) has not been 
varied because it is already known to be most probably small. A c omparison of the kno wn 
values of time lags in different Balmer lines of NGC 4151 allowed Shevchenkol f 1985a ) to 
conclude that the upper bound of the cavity radius Rq is 4-5 times less than the effective 
BL R radius. The ded uctio n that the cent r al cav ity is small is in agreement with conclusions 
by iMaoz et al.l (Il99ll ) and IXue fc ChengI (119981 ) . Note that the radius of the accretion disc 



in the centre of the nuclea r region is estimated to be equal to 0.6-2 It-days (iLyutyi 12005 



Sergeev et al.l . l2005l . l2006l ). Taking into account the uncertainty of this estimate and the 



uncertainty of the i?o estimates, the probable existence of this component of the nuclear 
region does not at all contradict the conclusions on the small relative size of the central 
cavity in the BLR. 

For the continuum flux we take the flux at the wavelengths of 4560-4640A ( the "4600A 
region), because this region corresponds to the shortest wavelengths at which iKaspi et al. 



(Il996h measured the continuum flux. The series of the observed values of the continuum 
flux are recalculated for presentation on the uniform time grid by means of cubic spline 
interpolation. 

The constant contribution to the integrated flux in Ha and H/3 due to the narrow com- 
ponents of the lines has been taken into account in the following way. We set Fa{t) = 
F^°'{t) + Fan and Fj3{t) = Fc'^it) + Fp^. The contributions of the narrow components Fan 
and Fpn are connected to each other via the constant Balmer decrement Dn = Fan/F^n- 
Therefore in the course of searching for the best model it is enough to vary the value of Fan] 
the value of is determined via the Balmer decrement. For each of the considered cases 
of orientation of the planes of clouds we have accomplished the modelling twice, namely, 
for the two reported va lues of the Balmer dec r ement in the narrow components: Dn = 4.47 
as given in Table 1 i n (IFerland fc Mushotzkyl . 119821 ) and Dn = 7.55 as given in Table 1 in 
(jSergeev et al.l . 1200 ll ) . The difference in the observed values of the decrement Dn may reflect 
either the difficulty in its evaluation or its probable long-term variability. 

The contribution of the stellar component to the observed continuum flux has been 
ta ken into accouii t by it s sub traction from the o bserved flux prior to modelling. According 
to Peterson et al. f 19951 ) and Kaspi et al. ( 1996 ). the stellar component contribution at the 
wavelength of 4600A and the aperture used at their observations is approximately equal to 
2.2 X 10"^^ erg cm'^ s"^ A'^. 

As it is known from observations (see, e.g., ICrenshaw et al.lll996l : iPeterson et al.ll2002l ). 



the light curves of an active galactic nucleus in the optical and UV continua can be rather dif- 
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ferent, and the relation between the continua is nonhnear. For the Seyfert galaxy NGC 5548, 
most studied in this respect, the slow components of variability in t he optical and UV con - 
tinua are connected by the power law Fopt oc F^y, where 7 ~ 0.56 (jPeterson et al.l . |2002| ). 
Basing on this relation, we find the real values of the s parameter from the values obtained 
in our modelling of the optical hght curves (we designate these values by Sopt) by means of 
the formula s = 7Sopt, where we set 7 = 0.6. 

The value of R does not depend on the line choice. The values of Sopt for different lines 
are generally different, and the same is true for a. The initial data for the iterations of the 
Levenberg-Marquardt algorithm have been taken randomly in the following limits: R — 
from 1 to 30 It-days; Sopt in the both lines — from 0.2 to 2.0; a — from 0.01 to 2.0; Fan — 
from to 30 X 10~^^ erg cm~^ s~^. 

In the case when the planes of clouds are orthogonal to the direction to the central 
source, the best-fit model light curves are presented in Fig. [Hlfor two different values of D^- 
The circles designate the observed values of the flux in optical continuura and the observed 
integrated emission-line fluxes according to the data in (IKaspi et al.lll996L Table 2). In both 
parts of the Figure, the continuous curve in the upper plot is the spline interpolation of the 
continuum flux. The flux is given in the units of 10~^^ erg cm~^ s~^ A~^. In the lower two 
plots, the best model light curves in Ha and H/5 are presented as continuous curves; here 
the integrated emission-line flux is given in the units of 10~^^ erg cm~^ s~^. 

In Fig. [71 the best model emission-line light curves are presented for the case when the 
planes of clouds are oriented randomly. 

To demonstrate the effect of variation of different parameters in determining the best-flt 
models, 'X^ /N (where N is the number of degrees of freedom) is shown in Figs. [HI [Ul and [TUl in 
dependence on the main parameters. Fig. [S] gives the dependence of the reduced on the 
BLR radius R with all other parameters assigned to their best-flt values. Figs. [91 and [TOl show 
the reduced in dependence on the s parameters for Ha and H/? with all other parameters 
assigned to their best-fit values. (Note that these graphs are presented here for illustrative 
purposes solely; for quantitative analysis, when offsetting one parameter, one should re- 
optimize all the other parameters while holding the one parameter fixed at the offset value.) 
In the presented graphs one can see that the dependence on R is characterized by a single 
well-defined minimum. The dependences of x^ on the s parameters have four minima each. 
The Levenberg-Marquardt algorithm finds the deepest one. 

In Tables [1] and [2], the obtained values of the model parameters, corresponding to the 
model emission-line light curves in Fig. [HI and Fig. [71 are presented together with the reduced 
X^ values of the models. In total, the results of our modelling give the following values of the 
radius of the BLR of NGC 4151: R= 11-14 It-days, with the uncertainty of 5-8 It-days. At 
Dn = 4.47, the best-fit values of the s parameter are: s ~ 0.6 for Ha and s ~ 0.85 for H/5. 
From comparison of Tables [T] and [2] one can see that at increasing decrement from 4.47 
up to 7.55 the difference in the computed values of the s parameter for Ha and B.(3 becomes 
les s, though does not seem to d isappear completely. 



Mushotzky &: Ferlandl (119841 ) elaborated photoionization models of stationary AGN opti- 



cal spectra. In the framework of the photoionization modelling they carried out calculations, 
in our equivalent terms representing the calculations of the functions L\{Fi) for the BLR 
clouds. According to the results of their calculations, s ~ 0.6 for Ha, and s ~ 0.8 for H/?. 
Thus, there exists a satisfact ory agreement of the results of our modelling with the data of 



Mushotzky fc Ferland (Il984[ ). especially in the case of = 4.47. Let us remark that the 



values of s can be nonconstant inside the BLR, varying from cloud to cloud, because they 
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Figure 6: The best model light curves of NGC 4151 in Ha and H/3 in the case when the planes 
of clouds are orthogonal to the direction to the central source. D^^ = 4.47 (a), — 7.55 
(b). The continuum flux is given in the units of 10^^"^ erg cm~^ s" -1 A-^ The emission-line 
fluxes are given in the units of 10~^^ erg cm~^ s~^. 
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Figure 7: The same as in Fig. [6], but the planes of clouds are oriented randomly. 
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Figure 8: The dependence of the reduced the BLR radius R with all other parameters 

assigned to their best-fit values. 



Table 1: The recovered values of the BLR parameters in the case of = 4.47 

Orthog. orient. Random orient. 

(Fig. Ek) (Fig. Uh) 

R 11 ±5 13 ±7 

s (sopt) 0.61 ± 0.08 (1.02 ± 0.13) 0.60 ± 0.08 (1.00 ± 0.14) 

Ha I F„n 19.76 ±0.97 19.55 ±1.14 

a 0.53 ±0.08 0.11 ±0.02 

s (sopt) 0.85 ± 0.05 (1.41 ± 0.09) 0.84 ± 0.06 (1.40 ± 0.10) 

H/3 I F/3n 4.42 ±0.22 4.37 ±0.26 

a 0.16 ±0.04 0.040 ±0.009 
X^N 2^97 2^98 

Fan and Fp^ are in the units of 10^^^ erg cm^^ s~^. 
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(b) s [Ha] 

Figure 9: The dependence of the reduced on the values of the s parameters for Ha and 
H/3 with all other parameters assigned to their best-fit values. The case when the planes of 
clouds are orthogonal to the direction to the central source. Dn = 4.47 (a), Dn = 7.55 (b). 
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Figure 10: The same as in Fig. [9l but the planes of clouds are oriented randomly. 
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Table 2: The same in the case of = 7.55 



Orthog. orient. Random orient. 

(Fig. Eb) (Fig. Cb) 



R 14 ±8 11 ±5 

s (sopt) 0.68 ±0.07 (1.13 ±0.12) 0.61 ± 0.10 (1.02 ± 0.16) 

Ha I F„n 20.96 ±0.86 20.22 ± 1.31 

a 0.56 ±0.13 0.11 ±0.01 

s (sopt) 0.71 ±0.05 (1.18 ±0.09) 0.66 ± 0.06 (1.10 ± 0.10) 

H/5 I Ff3n 2.78 ±0.11 2.68 ±0.17 

a 0.34 ±0.09 0.063 ± 0.009 

X^N 3.76 3.71 



depend on physical characteristics of the clouds. As a result of our modelling we obtain 
some "effective" values of s. 

The model of a homogeneous distribution of clouds implies that the covering factor is 
close to one, but, as it has been noted above, an interpretation of the adopted model is 
possible as a model with an exponential decrease of the cloud concentration with distance 
away from the centre; then the covering factor can be small. 



6 Discussion 



There are several ways of estimating the BLR size. Besides the reverberation and cross- 
correlation methods, discussed above, there exists a technique based on estimating the ion- 
ization pa rameter from modelling the s tationary emission-line spectra of AGN. Using this 
technique, iMushotzky fc Ferlandl (1l984l ) obtained an estimate of the radi us of the BLR of 



NGC 4151, equal to approximately 16 It-days. By a similar argument, ICassidy &: Raine 



(119971 ) found the inner and outer radii equal to 6 and 40 It-days in their theoretical model 
of the BLR of this galaxy . At t he same time when the photoionization estimate was made 
by IMushotzky fc Ferlandl f|l984), the reverberation estimate i? ~ 15 It-days was obtained 
independently by Shevchenkol ( 1984 ). in agreement with the photoionization estimate by 
Mushotzky fc Ferland ( 19841 ). This reverberation estimation was performed within the 

framework of the model of a homogeneous isotropic distribution of lir ie-emitting matter 

aroun d the central ionizing source, on t he ba sis of the observational data of lLyutyi fc Cherepashchuk 
(1l97ll ) and ICherepashchuk fc Lyutyil (1l973l ) on the time lags in the Ha line variations. 

Cross-correlation es timates are usual l y less than the "photoionizat ion" values. C ross- 
correlation analysis by Peterson fc Cota ( 1988 ) (see also discussion by Peterson 1988h. ac- 
compl ished on the basis of their own observational data and the data of lAntonucci Sz Cohen 
(119831 ) on variability in the lines H/3 and Hell A4686, gave ~ 6 It-days as the estimate for the 
radius of the BLR of NGC 4151. Similar cross-correlation estimates of the BLR size were 



recovered by lClavel et al.l (119901 ) on the basis of the lUE {International Ultraviolet Explorer) 
data on variability of the major UV lines: i? = 4 ± 3 It-days. These val ues correspond to the 
peak CCF time lags; the centroid ones are greater by about t wo days. IWandel et al.l (119991 ) 
find similar centroid CCF time lags, 4 ± 3 d, for the H/5 line. IClavel et al.l (Il990l ) note that 
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their cross-correlation estimates of the BLR size for NGC 4151 are an order of magnitude 
less than the typical "photoionization" estimates for Seyfert galaxies. 

By means of cross-correlation analysis of their own data, iKaspi et al.l (119961 ) found that 
the time lag of variations in the Ha and B.f3 lines in relation to continuum is 0-3 d; thus 
the cross-correlation estimate of the B LR radius is 0-3 It-day s. A ccording to t he mo dern 



analysis of these data accomplished by iMetzroth et al. 



mm andlBentz et all (120061 ). the 



Kaspi et al. (Il996h has no clear-cut 



value of the cross-correlation time lag for the data of 
statistical bounds. 

In total, the cross-correlation estimates of the BLR radius of NGC 4151 are all in the range 
of 0-6 It-days. Direct reverberation modelling, in comparison with the cross-correlation anal- 
ysis, give very different val ues of R sirn il ar to the given above "photoionization" estimates. 
According to the results of iMaoz et al.l (Il99ll ). who carried out reverberation modelling of 
light curves of NGC 4151 in Ha and H/5, the weighted-mean (by the local emission-line lumi- 
nosity of clouds) BLR radius ~ 16-18 It-days for the best found model, and the central cavity 
radius -Rq ~ 2 It-days. The linear character of the L](Fj) dependen ce was assumed, as in 
practically all modern research on this subject. Xue fc Chend (1998) nuni erical ly recovered 



the B LR transfer functions on the basis of the data of iMaoz et al.l (1199 ll ) and iKaspi et al. 
(119961 ). They obtained the following estimates: R ^ 10 It-days, Rq < 1 It-day. As men - 
tioned above, the reverberation estimate i? ~ 15 It-days was obtained in (jShevchenkd . 1 19841 ) . 
All these reverberation estimates are in agreement with our reverberation modelling results 
presented in Tables [H and El 

So, the known reverberation estimates of the BLR size of NGC 4151 are in agreement 
with "photoionization" estimates, and they all are much greater than the cross-correlation 
estimates. The strong difference between the BLR radii found by reverberation modelling, on 
one side, and its estimates following from cross-correlation analysis, on the other side, (10- 
18 It-days versus 0-6 It-days) underlines the conditional character of the cross-correlation 
estimates. Such a difference is no surprise: the size identified as the value of the observed 
time lag can be much (an order of magnitude) less than the true size of the BLR in It-days 
(Section H]). For example, if the cloud aggregate is uniform, the time lag of variation of a line 
with s ~ 1 with respect to an ionizing flare is small compared to the BLR radius R in light 
travel time units, and depends on R only weakly. The ultimate cause of this phenomenon is 
the degeneracy of relation ([5]) at s > 1. This degeneracy means that in practice there are no 
rigourous theoretical grounds to believe that the Aipeak value is mostly determined by the 
BLR size, if Atpeak is calculated for a typic al line (i.e., a line w ith s ~ 1). 

However, cross-correlation analysis by iKaspi et al.l (119961 ) of their observational data 
indicated that the cross-correlation time lag was small not only for H/? (the line with s ~ 1 
presumably), but for Ha as well (the line with s definitely less than one). To clarify this 
point, we have examined cross-correlations between the splined curve in continuum and our 
theoretical model light curves in Ha (presented in Figs. [HI and [71 linear trends have been 
subtracted prior to the analysis). The analysis has shown that the cr oss-correlatio r i time lag 
is 0-3 d, practically equal to the cross-correlation time lag found by iKaspi et al.l (1l996l ) for 
the observed emission-line light curves; the computed cross-correlation functions are shown 
in Fig. [m So, substitution of the observed emission-line light curve by a theoretical one 
calculated for a definite value of R (equal to 11-14 It-days) in the transfer function does 
not change the cross- correlation estimate of R (0-3 It-days). It is much less than the value 
adopted as the parameter of the transfer function. This means that the cross-correlation 
method may provide inadequate results not only at s ~ 1, but at s < 1 as well. The cause 
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Figure 11: The computed cross-correlations between the sphned curve in continuum and our 
theoretical model light curves. The case when the planes of clouds are orthogonal to the 
direction to the central source (a) and the case when they are oriented randomly (b). 
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is not all-together clear, but one may speculate that it is related to the strong dependence 
of Atpcak on the variability time scale (Fig. E]) and/or to the small amplitudes of variability 
in this particular set of observational data. Let us underline that, contrary to the cross- 
correlation analysis, the reverberation modelling turns out to be immune, as we have seen, 
to these unfavourable conditions, and provides adequate values of R. 

A hypothes i s on possible variability of t he BLR size of NGC 4151 was put forward 
by iKaspi et al.l (119961 ) , iPeterson et al.l (120021 ) and other researchers on the basis of cross- 
correlation analysis of optical spectral variability data at different time interval s of ob serva- 
tions. Our reverberation modelling of the Ha light curve data of iKaspi et al.l (119961 ) gives 
the value of the BLR radius matching the majority of the BLR size estimates of other au- 
thors. This removes necessity in any special physical interpretation of the small value of the 
cross-correla t ion ti me lag in Ha for these light curve data. In particular, the hypothesis by 



Kaspi et al.l (119961 ) that the physical size of the BLR at the moment of their observations 



was an order of magnitude less than usually is not necessitated. 



7 Conclusions 



We have studied how the nonlinearity in the "Li-Fi" relation (the emission-line luminosity, 
Li, of the BLR cloud in dependence on the ionizing continuum flux, Fi, incident on the 
cloud) can be taken into account in estimating the size of the BLR in active galactic nuclei 
by means of "reverberation" methods. We have shown that the BLR size estimates obtained 
by cross-correlation peaks of emission-line and continuum light curves can be much (up 
to an order of magnitude) less than the values obtained by reverberation modelling. This 
has been demonstrated by means of abstract representative numerical cross-correlation and 
reverberation experiments with model continuum flares and emission-line transfer functions 
and by means of practical reverberation modelling of the observed emission-line variability 
of NGC 4151. The modelling of the observed light curves of N GC 4151 in Ha ari d H/3 has 
been accomplished on the basis of th e observatio n al da t a bv iKaspi et al.l (119961 ) and the 



theoretical broad- line region model by IShevchenkd (119841 . Il985al ). 

In the abstract representative numerical cross-correlation and reverberation experiments 
with model continuum flares and emission-line transfer functions, we have found that the 
value of the cross-correlation peak time lag Atpeak for s > 1 is small in comparison with the 
BLR size R expressed in the light-travel time units and depends on R only weakly. We have 
shown that in the case of s > 1 the effect of the ionizing flare duration on the Atpeak value is 
far greater than that of the BLR radius. In other words, the BLR radius has little effect on 
the measured value of the Atpeak value in the mathematically degenerate but observationally 
most common case of s = 1; the role of the timescale of variability is far greater. Therefore, 
the lines with s ~ 1 seem to be of little help in determining the size of the BLR by means 
of estimating the cross-correlation peak time lag. 



T he presence of a noticeable t ime la g of variations of NGC 4151 in Ha ( ILyutyi &: Cherepashchukl . 
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Cherepashchuk fc Lyutvil. 19731 ) and significantly shorter time lags in other Balmer 



lines (lAntonucci fc Cohenl . Il983l ) with respect to variations i r i the optical continuum has 



been attributed, in agreement with conclusions by IShevchenkd (119841 . Il985al ). to the effect of 
essential nonlinearity in the "Li-F;" relation for Ha. The low value of the power-law index, 
s ~ 0.6, distinguishes this line from the other Balmer lines. 

The values of the model parameters of the BLR of NGC 4151 have been estimated. In 
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particular, estimates of the BLR radius have been made. Our r everberation modelhng of 



the emission-hne variabihty based on the observational data by iKaspi et al.l (Il996l ) gives 
values of the BLR radius agreeing with the majority of its known "reverberation" and "pho- 
toionization" estimates. Much smaller R values obtained by means of the cros s-correlation 



metho d have been shown to be an artifact of this method. The hypothesis by iKaspi et al. 



(119961 ) that the size of the BLR of NGC 4151 at the time interval of their observations were 
an order of magnitude less than usually is not necessitated. 

Concluding, a power-law emission-line response model and simple spherically symmet- 
ric thick geometries of the BLR cloud distribution, taken here as a basis for modelling the 
emission- line AGN variability, gives the size of the BLR of NGC 4151 equal to 11-14 It-days, 
with the uncertainty of 5-8 It-days. This agrees satisfactorily with BLR size estimates in 
photoionization models fitting emission line strengths in the mean spectrum of NGC 4151. 
Much shorter time lags found in the cross-correlation analysis of the emission-line and con- 
tinuum light curves of NGC 4151 correspond to the size of a smaller emission-line region 
that is reverberating. 
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